Partition satellite mission planning method for earliest completion of regional target coverage

ABSTRACT

The present invention discloses a partition satellite mission planning method for earliest completion of regional target coverage, which includes the following steps: 1. partitioning a to-be-observed rectangular region; 2. allocating observation resources to different regions and selecting certain coverage opportunities and coverage modes thereof to minimize completion time corresponding to a formed coverage solution. The present invention can obtain a satisfactory solution capable of finishing complete coverage of regional targets as early as possible with sufficient observation resources through suitable calculation resources and achieve balance between the consumption of calculation resources and optimality of solutions, so that satellite resources can be fully used to conduct the quick and effective coverage search for important regional targets in a real environment.

CROSS-REFERENCE TO RELATED APPLICATIONS

This application claims the benefit of priority from Chinese Patent Application No. 201910957601.4, filed on Oct. 10, 2019. The content of the aforementioned applications, including any intervening amendments thereto, is incorporated herein by reference in its entirety.

TECHNICAL FIELD

The present invention belongs to the field of regional target mission planning for earth observation satellites, and particularly relates to a partition satellite mission planning method for earliest completion of regional target coverage.

BACKGROUND OF THE PRESENT INVENTION

Earth observation satellites (EOS) or imaging satellites (abbreviated as satellites hereinafter) have a main function of observing land, oceans and atmosphere through satellite borne sensors (such as visible light cameras and multi-spectral cameras). Observation requirements are proposed by users from various fields and departments and collected by a ground control center of the satellites. The control center comprehensively formulates an imaging coverage plan of each satellite according to the observation requirements in combination with the application conditions of satellite resources, and generates measurement and control instructions to be transmitted to the satellites via a ground measurement and control station. The satellites make corresponding actions after receiving the instructions to image a specified region. The image data is temporarily saved in a satellite borne hard drive. When the satellites can communicate with the ground station, the image data is transmitted to the ground station. In this process, the link that the ground control center formulates the satellite imaging plan is called satellite mission planning, which is one of key links in the entire satellite application management process.

In the traditional satellite application modes, the satellites formulate the plans separately and perform the imaging missions separately. No collaborative observation is carried out among different satellites. As the number of satellites increases, it is possible to use multiple satellites to collaboratively perform the imaging observation for earth. The collaborative regional imaging of multiple satellites is actually needed.

The satellite can only photograph a strip region of limited length and width at one pass. If the to-be-observed region is large, it is difficult for the satellite to completely observe the entire region at one pass. If the user urgently needs the image data of the region, multiple satellites can be used to perform the collaborative imaging at multiple passes. Since tracks of sub-satellite points of each satellite are not necessarily parallel, the imaging strips may be overlapped. If the satellites are not properly arranged, a large number of strips may be overlapped, causing repeated observation of some regions, while some regions are not observed, which is not conducive to the rational utilization of resources. In order to better use the existing satellite coverage resources, a rational plan should be formulated to solve the above contradiction. It is a typical optimization problem of operations research, which needs to be urgently solved. When the to-be-observed region with this problem is large, a conventional regional target planning method may consume a great amount of calculation resources due to the large scale and high complexity of the problem. This is unfavorable for the application of the method in actual engineering environments.

SUMMARY OF THE PRESENT INVENTION

In order to solve the defects in the prior art, the present invention provides a partition satellite mission planning method for earliest completion of regional target coverage so as to obtain a satisfactory solution capable of completing the complete coverage of regional targets as early as possible with sufficient observation resources through suitable calculation resources, and to achieve balance between the consumption of calculation resources and the optimality of solutions, so that the satellite resources can be fully used to conduct quick and effective coverage search for important regional targets in a real environment.

To solve the technical problems, the present invention adopts the following technical solutions:

Compared with the prior art, the present invention has the following beneficial effects:

1. According to a “divide-and-rule” concept, the present invention proposes a solving strategy based on partition, that is, a large to-be-observed region is divided into a plurality of small subregions; coverage resources are allocated to each subregion; each subregion is solved separately; and the solutions of all subregions are combined to form overall coverage solution, thereby reducing the problem scale and alleviating the optimization difficulty. Different allocation solutions of coverage resources may lead to different overall coverage solutions, so that a simulated annealing method is proposed to search an optimal resource allocation solution.

2. A grid discretization method is used to process the to-be-observed subregions, and the original problem is converted to a problem of covering each grid, so that the coupling with the calculation geometry can be reduced. At the same time, based on grids, a longest basic coverage mode set of each observation opportunity is built as an optional coverage mode set, and decision variables are changed from values in a continuous space to values in a discrete space, thereby reducing the calculation amount.

3. The present invention utilizes a two-phase heuristic algorithm to find a desirable coverage solution for each discretized subregion. At a first phase, a dynamic greed rule is used to find a desirable feasible coverage solution capable of being quickly obtained under the current environment. At a second phase, after a coverage opportunity with maximum end time in the selected coverage opportunities is removed on the basis of the existing solutions, the coverage solution is re-planned to improve the quality of the solutions. The two-phase heuristic algorithm is a polynomial time algorithm. When the total number of the coverage modes is in a rational range, the algorithm is high in solving speed, easy to implement, high in interpretability, stable and reliable, so that in the simulated annealing process of the method, a desirable coverage solution can be rapidly and effectively found to increase the search speed for a solution space of allocation solutions.

BRIEF DESCRIPTION OF THE DRAWINGS

FIG. 1 is a flow chart of a partition satellite mission planning method for earliest completion of regional target coverage of the present invention;

FIG. 2a is a schematic diagram of left lower corner oblique coverage of a region target of the present invention;

FIG. 2b is a schematic diagram of right lower corner oblique coverage of the region target of the present invention;

FIG. 2c is a schematic diagram of multi-satellite collaborative coverage of the present invention;

FIG. 3a is a schematic diagram of formation of coverage modes of the present invention;

FIG. 3b is a schematic diagram of a longest coverage mode length of the present invention;

FIG. 3c is a schematic diagram of side-sway of satellites of the present invention;

FIG. 3d is a schematic diagram of a maximum roll angle of the present invention;

FIG. 4a is a schematic diagram of a longest basic coverage mode of a right lower corner oblique type of the present invention;

FIG. 4b is a schematic diagram of the longest basic coverage mode of a left lower corner oblique type of the present invention;

FIG. 5a is a schematic diagram of calculation of a roll angle under one condition of the present invention;

FIG. 5b is a schematic diagram of calculation of a roll angle under another condition of the present invention;

FIG. 6a is a schematic diagram of a feasible upper grid condition 1 of the present invention;

FIG. 6b is a schematic diagram of a feasible upper grid condition 2 of the present invention;

FIG. 6c is a schematic diagram of a feasible upper grid condition 3 of the present invention;

FIG. 6d is a schematic diagram of a feasible upper grid condition 4 of the present invention.

DETAILED DESCRIPTION OF THE PRESENT INVENTION

As shown in FIG. 2a , FIG. 2b and FIG. 2c , in the present embodiment, a partition satellite mission planning method for earliest completion of regional target coverage is used to select a subset from sufficient coverage opportunity set S and to plan corresponding coverage modes for completely covering a to-be-observed region R. The optimization purpose of this method is to minimize the completion time. The coverage opportunity refers to an opportunity for a satellite to perform the coverage observation when passing over the to-be-observed region. As shown in FIG. 2c , an area of the to-be-observed region R is large. The single coverage opportunity can only cover a small part, so that a plurality of coverage opportunities are needed for collaborative observation.

As shown in FIG. 3a , FIG. 3b , FIG. 3c and FIG. 3d , each coverage opportunity s∈s has own attribute and includes sub-satellite point track o_(s), height h_(s) of the satellite to the ground and maximum power-on time m_(s), corresponding to a maximum length of coverage mode r_(s), maximum roll angle v_(s) and camera field angle w_(s). The satellite orbits around the Earth, and a track formed by sub-satellite points on the ground should be a curve actually, but it can be approximated to a straight line in a small range inside the to-be-observed region. In order to simplify the problem, the sub-satellite point track o_(s) of the coverage opportunity is considered as the straight line in the present embodiment of the present invention. As shown in FIG. 2a and FIG. 2b , each coverage opportunity has a start time α(s) and an end time β(s) for the to-be-observed region R.

The so-called coverage mode c is a rectangular strip photographed by the satellites and determined by a certain roll angleand on/off time selected by one coverage opportunity.

A coordinate system o-xy is established by adopting any vertex of the to-be-observed rectangular region R as an origin o and two sides adjacent to the origin as axis x and axis y. The coordinates of four vertexes of the to-be-observed rectangular region R are marked as left upper corner vertex LU(R), left lower corner vertex LD(R), right upper corner vertex RU(R) and right lower corner vertex RD(R) respectively. The partition satellite mission planning method includes the following steps:

Step 1: a to-be-observed rectangular region R is partitioned.

Step 1.1: length L_(P) and width W_(P) of each subregion are set.

Step 1.2: the length L_(R) and width W_(R) of the to-be-observed rectangular region R are calculated according to the coordinates of four vertexes of the to-be-observed rectangular region R.

Step 1.3: L_(P) is divided by L_(R) to obtain length section number D_(L) and remainder R_(L) of the subregion. If the remainder

${R_{L} \geq {\frac{D_{L}}{{2D_{L}} + 1}L_{p}}},{D_{L} + 1}$

is assigned to D_(L). W_(p) is divided by W_(R) to obtain width section number D_(W) and remainder R_(W) of the subregion, and if the remainder

${R_{W} \geq {\frac{D_{W}}{{2D_{W}} + 1}W_{p}}},{D_{W} + 1}$

is assigned to D_(W).

Step 1.4: L_(R) is divided by D_(L) to obtain length L_(f) of the final subregion, and W_(R) is divided by D_(W) to obtain width W_(f) of the final subregion.

Step 1.5: according to the length L_(f) and width W_(f) of the final subregion, the to-be-observed rectangular region R is equally partitioned, thereby obtaining the coordinates of four vertexes of each subregion.

Step 1.5.1: let i=1, the left upper corner coordinate LU(R) of the to-be-observed rectangular region R is used as the left upper corner coordinate LU(n_(i)) of the i^(th) subregion n_(i).

Step 1.5.2: the abscissa value of the left upper corner coordinate LU(n_(i)) of the i^(th) subregion n_(i) is added with L_(f) to obtain the right upper corner coordinate RU(n_(i)) thereof. The ordinate value of the left upper corner coordinate LU(n_(i)) is added with W_(f) to obtain the left lower corner coordinate LD(n_(i)) thereof. The ordinate value of the right upper corner coordinate RU(n_(i)) is added with f to obtain the right lower coordinate RD(n_(i)).

Step 1.5.3: if the right lower corner coordinate RD(n_(i)) of the i^(th) subregion n_(i) is equal to the right lower corner coordinate RD(R) of the to-be-observed rectangular region R, step 1.5.6 is carried out; otherwise, step 1.5.4 is carried out.

Step 1.5.4: if the abscissa value of the right upper corner coordinate RU(n_(i)) of the i^(th) subregion n_(i) is equal to the abscissa value of the right upper corner coordinate RU(R) of the to-be-observed rectangular region R, LD(n_((i−D) _(L) ₊₁₎) is used as the left upper corner coordinate LU(n_(i+1)) of the (i+1)^(th) subregion n_(i+1); otherwise, the right upper corner coordinate RU(n_(i)) of the i^(th) subregion n_(i) is used as the left upper corner coordinate LU(n_(i+1)) of the (i+1)^(th) subregion n_(i+1).

Step 1.5.5: after i+1 is assigned to i, step 1.5.2 is carried out to construct the coordinates of four vertexes of a next subregion.

Step 1.5.6: the partitioning is ended, and the coordinates of four vertexes of each subregion are outputted.

Step 1.6: a set N of the subregions is obtained, and the total number of the subregions is |N|=D_(L)×D_(W).

Step 2: observation resources are allocated to different regions, and a certain subset of coverage opportunities and corresponding coverage modes are selected to minimize the completion time of a formed coverage solution. That is, different coverage opportunities are responsible for covering different subregions, while different allocation solutions may finally lead to different coverage solutions which correspond to different completion time. Therefore, a simulated annealing method is used to search an optimum allocation solution here.

Step 2.1: the coverage opportunities are randomly allocated to each subregion to form an initial current allocation solution Y={y^(s)|∀s∈S} wherein the value of y_(s) represents the subregion where the coverage opportunity s is allocated.

Step 2.2: the coverage opportunity set S*(Y) selected under the current allocation solution Y, a desirable coverage mode set C*(Y) and corresponding completion time F(Y) are calculated.

Step 2.2.1: different subregions and the set of coverage opportunities allocated to different subregions are used as inputs to construct grids, and a two-phase heuristic algorithm is used to calculate the coverage opportunity set selected for different subregions as well as the desirable coverage solution and corresponding completion time.

Step 2.2.1.1: the i^(th) subregion is divided into a plurality of consistent square grids to obtain a grid set J_(i), wherein a side length of each grid should be minimized. The existing grids are ensured to be completely covered by the coverage modes.

If the to-be-observed region was not partitioned into subregions, the whole large region may form a large number of grids A large number of coverage modes may be generated based on these grids, which consumes vast calculation resources. Some large-scale calculation examples cannot even obtain a solution within a rational time, so that embodiments of the present invention are proposed on the basis of solving strategies based on the partition.

Step 2.2.1.2: according to the divided grid set J_(i), a set of longest basic coverage modes is generated for each coverage opportunity in set S_(i) allocated to the i^(th) subregion n_(i), and jointly form a main set C_(i) of the optional coverage modes.

As shown in FIG. 4a and FIG. 4b , the longest basic coverage mode refers to that the length is equal to the maximum length r_(s), and there are two grids that satisfy the following two conditions:

{circle around (1)} Both grids are completely covered by the coverage mode, that is, the four vertexes are within the area covered by the coverage mode;

{circle around (2)} one grid has a vertex located on a left border of the coverage mode and is named left grid, and the other grid has a vertex located on an upper border of the coverage mode and is named upper grid. The left grid and the upper grid can be the same.

Each coverage opportunity can select infinitely many on/off time and roll angles, so each coverage opportunity may actually correspond to an infinite number of coverage modes. However, all coverage modes cannot be enumerated. Therefore, the solving method of the problem is difficult to design on the basis of the universal set of all coverage modes. The longest basic coverage mode is a special coverage mode and can substitute a majority of non-longest basic coverage modes and is limited in quantity. A majority of non-longest basic coverage modes may be slightly adjusted to be the corresponding longest basic coverage modes without reducing the number of completely covered grids.

Step 2.2.1.2.1: if S_(i)=Ø, step 2.2.1.2.4 is carried out; otherwise, step 2.2.1.2.2 is carried out.

Step 2.2.1.2.2: one of coverage opportunities s∈S_(i) is selected, and according to the type of the coverage opportunity s, methods in step 2.2.1.2.2.1 to step 2.2.1.2.2.10 are used to construct a longest basic coverage mode set C_(s) based on the grid set J_(i) for the coverage opportunity s.

Step 2.2.1.2.2.1: the type of the coverage opportunity s is judged according to the sub-satellite point track o_(s) of the coverage opportunity s, and the types include a right lower corner oblique type and a left lower corner oblique type.

Step 2.2.1.2.2.2: let set C_(s)=Ø.

Step 2.2.1.2.2.3: subsets J_(s) are selected from the grid set J_(i), and each grid in the subsets J_(s) may be used as the left grid of a longest basic coverage mode of the coverage opportunity s.

Step 2.2.1.2.2.3.1: let J_(s)=Ø.

Step 2.2.1.2.2.3.2: the set J_(i) is traversed, and the camera roll angle

${\sigma_{s}(u)} = {{\arctan \frac{d_{s}^{u}}{h_{s}}} - \frac{w_{s}}{2}}$

of the coverage opportunity s is calculated when each grid u is used as the left grid of the coverage mode. If the coverage opportunity is the right lower corner oblique type, as shown in FIG. 4a , d_(s) ^(u) is a distance between the left lower corner vertex LD(u) of the grid u and the sub-satellite point track o_(s) of the coverage opportunity s. If the coverage opportunity is the left lower corner oblique type, as shown in FIG. 4b , d_(s) ^(u) is a distance between the left upper corner vertex LU(u) of the grid u and the sub-satellite point track o_(s) of the coverage opportunity s. If |σ_(s)(u)|≤v_(s), the grid u is added into the subset J_(s) of the left grid.

Step 2.2.1.2.2.4: if J_(s)=Ø, step 2.2.1.2.2.10 is carried out; otherwise, after one grid p in J_(s) is selected, step 2.2.1.2.2.5 is carried out.

Step 2.2.1.2.2.5: according to the grid p, the grid subset J_(s) ^(p) is selected from the grid set J_(i), and when p is used as the left grid, each grid in J_(s) ^(p) may be used as the upper grid to form a longest basic coverage mode of the coverage opportunity s together with p.

Step 2.2.1.2.2.5.1: let J_(s) ^(p)=Ø.

Step 2.2.1.2.2.5.2: a strip width η_(s)(p) of the coverage mode when the grid p is used as the left grid is calculated by using the formula (1):

$\begin{matrix} {{\eta_{s}(p)} = \left\{ \begin{matrix} {{d_{s}^{p} - {h_{s} \cdot {\tan \left( {{\arctan \frac{d_{s}^{p}}{h_{s}}} - w_{s}} \right)}}},{{\arctan \; \frac{d_{s}^{p}}{h_{s}}} > w_{s}}} \\ {{d_{s}^{p} + {h_{s} \cdot {\tan \left( {w_{s} - {\arctan \frac{d_{s}^{p}}{h_{s}}}} \right)}}},{{\arctan \; \frac{d_{s}^{p}}{h_{s}}} \leq w_{s}}} \end{matrix} \right.} & (1) \end{matrix}$

As shown in FIG. 5a and FIG. 5b , when the grid u is used as the left grid of the coverage mode, a calculation formula of the camera roll angle σ_(s)(u) of the coverage opportunity s and the strip width η_(s)(u) of the coverage mode can both be obtained by elementary geometry. OB is a vertical line, and the length of OB is the height of the satellite to the ground, i.e. h_(s). The length of AB is a distance between the straight lines L_(u) and o_(s), i.e. d_(s) ^(u). Apparently, ∠ABO is a right angle, so it is easy to get:

${{\angle \; {AOB}} = {\arctan \frac{d_{s}^{u}}{h_{s}}}},$

and the camera roll angle

${\sigma_{s}(u)} = {{\arctan \frac{d_{s}^{u}}{h_{s}}} - {\frac{w_{s}}{2}.}}$

In FIG. 5a , ∠AOC is a field angle, i.e. w_(s), while ∠AOB>w_(s). The width of the strip is the length of AC at this moment.

${AC} = {{d_{s}^{u} - {BC}} = {{d_{s}^{u} - {{h_{s} \cdot \tan}\; \angle \; {BOC}}} = {{d_{s}^{u} - {h_{s} \cdot {\tan \left( {{\angle \; {AOB}} - w_{s}} \right)}}} = {d_{s}^{u} - {h_{s} \cdot {\tan \left( {{\arctan \frac{d_{s}^{u}}{h_{s}}} - w_{s}} \right)}}}}}}$

In FIG. 5b , ∠AOD is a field angle, i.e. w_(s). OC is an angular bisector of ∠AOD, i.e.

${{\angle AOC} = \frac{w_{s}}{2}}.$

The width of the strip is length of AD at this moment.

${{\angle \; {AOB}} = {\arctan \; \frac{d_{s}^{u}}{h_{s}}}},{{\angle \; {DOB}} = {w_{s} - {\angle \; {AOB}}}},$

so it can obtain: BD=h_(s)·tan ∠DOB, and it can further obtain:

${AD} = {{d_{s}^{u} + {BD}} = {{d_{s}^{u} + {h_{s} \cdot {\tan \angle {BOD}}}} = {{d_{s}^{u} + {h_{s} \cdot {\tan\left( {w_{s} - {\angle {AOB}}} \right)}}} = {d_{s}^{u} + {h_{s} \cdot {\tan \left( {w_{s} - {\arctan \frac{d_{s}^{u}}{h_{s}}}} \right)}}}}}}$

Step 2.2.1.2.2.5.3: the set J_(i) is traversed. Whether one grid t satisfies the following four conditions is checked. If all conditions are satisfied, the grid t is added into J_(s) ^(p); otherwise, the next grid is checked.

If the coverage opportunity is the right lower corner oblique type, the four conditions are:

{circle around (1)} as shown in FIG. 6a , the straight line L_(p)∥o_(s) is made through the left lower corner vertex of the grid p, the straight line L_(t)′∥o_(s) is made through the right upper corner vertex of t, and the distance between the straight lines L_(p) and L_(t)′ cannot exceed the width η_(s)(p) of the strip.

{circle around (2)} as shown in FIG. 6b , the straight lines L_(p)′⊥o_(s) and L_(t)⊥o_(s) are made respectively through the left upper corner vertexes of the grids p and t, and L_(t) is above L_(p)′.

{circle around (3)} as shown in FIG. 6c , the straight line L _(t)∥o_(s) is made through the left lower corner vertex of the grid t, and the straight line L _(t) is at the right side of L_(p).

{circle around (4)} as shown in FIG. 6d , the straight line L _(p)⊥o_(s) is made through the right lower corner vertex of the grid p, and the distance between the straight lines L _(p) and L_(t) cannot exceed the length l_(s) of the strip. “∥” indicates a parallel relation, and “⊥” indicates a vertical relation.

If the coverage opportunity is the left lower corner oblique type, the four conditions are:

{circle around (1)} The straight line L_(p)∥o_(s) is made through the left upper corner vertex of the grid p, the straight line L_(t)′∥o_(s) is made through the right upper corner vertex of t, and the distance between the straight lines L_(p) and L_(t)′ cannot exceed the width η_(s)(p) of the strip.

{circle around (2)} The straight lines L_(p)′⊥o_(s) and L_(t)⊥o_(s) are made respectively through the right upper corner vertexes of the grids p and t, and L_(t) is above L_(p)′.

{circle around (3)} The straight line L _(t)∥o_(s) is made through the left upper corner vertex of the grid t, and the straight line L _(t) is at the right side of L_(p).

{circle around (4)} The straight line L _(p)⊥o_(s) is made through the left lower corner vertex of the grid p, and the distance between the straight lines L _(p) and L_(t) cannot exceed the length l_(s) of the strip. “∥” indicates a parallel relation, and “⊥” indicates a vertical relation.

Apparently, the grid p also satisfies the four conditions.

Step 2.2.1.2.2.6: if J_(s) ^(p)=Ø, step 2.2.1.2.2.9 is carried out; otherwise, after one grid q in J_(s) ^(p) is selected, step 2.2.1.2.2.7 is carried out.

Step 2.2.1.2.2.7: based on the grids P and q, the longest basic coverage c_(s)(p,q) belonging to the coverage opportunity s is constructed.

Step 2.2.1.2.2.7.1: if the coverage opportunity is the right lower corner oblique type, the straight line L_(p)∥o_(s) is made through the left lower corner vertex LD(p) of the grid p, and the straight line L_(q)⊥o_(s) is made through the left upper corner vertex LU(q) of the grid q. If the coverage opportunity is the left lower corner oblique type, the straight line L_(p)∥o_(s) is made through the left upper corner vertex LU(p) of the grid p, and the straight line L_(q) ⊥o_(s) is made through the right upper corner vertex RU(q) of the grid q.

Step 2.2.1.2.2.7.2: the strip width η_(s)(p) of the coverage mode when the grid p is used as the left grid is calculated.

Step 2.2.1.2.2.7.3: the straight line L _(p) is made so as to satisfy:

{circle around (1)} L _(p)∥L_(p);

{circle around (2)} (the distance between L _(p) and L_(p) is equal to η_(s)(p);

{circle around (3)} L _(p) is at the right side of L_(p).

Step 2.2.1.2.2.7.4: the straight line L _(q) is made so as to satisfy:

{circle around (1)} L_(q)∥L_(q);

{circle around (2)} (the distance between L _(q) and L_(q) is equal to l_(s);

{circle around (3)} L _(q) is at the lower side of L_(q).

Step 2.2.1.2.2.7.5: four straight lines of L_(p), L _(p), L_(q), L _(q) are used as four sides of the rectangle, and intersections of the four straight lines are used as four vertexes of the rectangle to construct the coverage mode c_(s)(p,q).

Step 2.2.1.2.2.8: c_(s) (p,q) is added into C_(s), the grid q is removed from J_(s) ^(p), and then step 2.2.1.2.2.6 is carried out.

Step 2.2.1.2.2.9: the grid p is removed from J_(s), and then step 2.2.1.2.2.4 is carried out.

Step 2.2.1.2.2.10: the coverage mode set C_(s) is outputted.

Step 2.2.1.2.3: the coverage mode opportunity s is removed from S_(i), and then step 2.2.1.2.1 is carried out.

Step 2.2.1.2.4: let, the coverage mode set C_(i) is outputted.

Step 2.2.1.3: by utilizing the grid set J_(i) separated from i_(th) subregion n_(i) as well as the coverage opportunity set S_(i)={s|y_(s)=n_(i), s∈S} allocated to i^(th) subregion n_(i) and optional coverage mode set C_(i) as inputs, the coverage opportunity subset S_(i)* selected by the i^(th) subregion n_(i) is determined by the two-phase heuristic algorithm, and the coverage mode of each coverage opportunity is selected, so that the desirable coverage mode set C_(i)* is obtained, and the completion time F_(i) thereof can be obtained at the same time.

Step 2.2.1.3.1: by utilizing the grid set J_(i), the coverage opportunity set S_(i) and the corresponding coverage mode set C_(i) as inputs, the selected coverage opportunity subsets S_(i)*, the selected coverage mode sets C_(i)* and the number A_(i) of corresponding coverable grids can be calculated with the heuristic algorithm.

Step 2.2.1.3.1.1: all coverage opportunities in S_(i) are sorted in a non-descending order according to the end time to obtain

$\overset{\rightharpoonup}{S_{i}} = {\left\{ {s_{1},s_{2},\ldots,s_{|\overset{\rightharpoonup}{s_{i}|}}}\  \right\}.}$

Step 2.2.1.3.1.2: the state of all grids in J_(i) is set as “uncovered”, let g=0, and C_(i)*=Ø.

Step 2.2.1.3.1.3: g+1 is assigned to g, the coverage opportunity s_(g) is extracted from {right arrow over (S_(i))}, and the optional coverage mode set C_(s) _(g) is traversed. For each coverage mode C_(s) _(g) , the grid capable of being completely covered by C_(s) _(g) and with the state of “uncovered” is selected from the grid set J_(i) and marked as valid grid sets V(c_(s) _(g) ), while the number is marked as |V(c_(s) _(g) )|.

Step 2.2.1.3.1.4: the coverage mode

$c_{s_{g}}^{*} = \left. {\arg \max\limits_{\forall{c_{s_{g}} \in C_{s_{g}}}}} \middle| {V\left( c_{s_{g}} \right)} \right|$

covering the most valid grids is selected from C_(s) _(g) and used as the coverage mode selected by the coverage opportunity s_(g). After c_(s) _(g) * is added into C_(i)*, the state of all grids in the valid grid sets V(c_(s) _(g) ) of the coverage mode c_(s) _(g) * is changed as “covered”.

Step 2.2.1.3.1.5: the number of all grids with the state of “covered” in J_(i) is calculated and marked as A_(i).

Step 2.2.1.3.1.6: if A_(i)=|J_(i)| or g=|{right arrow over (S_(i))}|, step 2.2.1.3.1.7 is carried out; otherwise, step 2.2.1.3.1.3 is carried out.

Step 2.2.1.3.1.7: the coverage opportunity subset S_(i)*={s₁, s₂, . . . , s_(g)} as well as the coverage mode set C_(i)* thereof and the number A_(i) of covered grids are outputted.

Step 2.2.1.3.2: if A_(i)<|J_(i)|, step 2.2.1.3.4 is carried out; otherwise, step 2.2.1.3.3 is carried out.

On the premise of sufficient coverage opportunities, the coverage opportunities that are randomly allocated to each subregion generally can completely cover the subregion, that is, A_(i)=|J_(i)|. If some subregions cannot be completely covered, they may be improved in neighbor allocation solutions constructed by an insert operator and an exchange operator at the simulated annealing phase.

Step 2.2.1.3.3: S_(i)* and C_(i)* are optimized by using a promotion algorithm.

Step 2.2.1.3.3.1: one coverage opportunity

$s_{m} = {\arg {\max\limits_{\forall{s \in S_{i}^{*}}}{\beta (s)}}}$

with maximum end time in S_(i)* is removed from S_(i)* to obtain the coverage opportunity set S_(i) .

Step 2.2.1.3.3.2: by utilizing the grid set J_(i), the coverage opportunity set S_(i) and the optional coverage mode set C_(i) ={c|∀c∀C_(s), ∀s∈S_(i) } corresponding to S_(i) as inputs, the selected coverage opportunity subsets S_(i)*, the selected coverage mode set C_(i)*, and the number A_(i) of the corresponding coverable grids can be calculated with the heuristic algorithm.

Step 2.2.1.3.3.2.1: the state of all grids in J_(i) is set as “uncovered”, and C_(i)*=Ø.

Step 2.2.1.3.3.2.2: let g=1.

Step 2.2.1.3.3.2.3: if g>|C_(i) |, step a.3.3.2.5 is carried out; otherwise, step a.3.3.2.4 is carried out.

Step 2.2.1.3.3.2.4: the g ^(th) coverage mode in C_(i) is selected, and the grids capable of being completely covered by the g ^(th) coverage mode and with the state of “uncovered” are selected from the grid set J_(i) and marked as valid grid set V(g), while its number is marked as |V(g)|. g+1 is assigned to g. Thereafter step a.3.3.2.3 is carried out.

Step 2.2.1.3.3.2.5: if C_(i) =Ø, step a.3.3.2.8 is carried out; otherwise, step a.3.3.2.6 is carried out.

Step 2.2.1.3.3.2.6: the coverage mode

${\overset{\_}{g}}^{*} = \left. {\arg \max\limits_{\forall{\overset{\_}{g} \in \overset{\_}{C_{i}}}}} \middle| {V\left( \overset{\_}{g} \right)} \right|$

with most valid grids in C_(i) is selected and used as the coverage mode selected by the corresponding coverage opportunity s*. g* is added into C_(i)*, and the state of all grids in the valid grid set V(g*) of the coverage mode g* is changed to be “covered”.

Step 2.2.1.3.3.2.7: the coverage mode set C _(s*) corresponding to the coverage opportunity s* is removed from the optional coverage mode set C_(i) , and then step a.3.3.2.2 is carried out.

Step 2.2.1.3.3.2.8: the coverage mode set C_(i)*, selected by each coverage opportunity is outputted, i.e. the final coverage solution. The number of all “covered” grids in J_(i) is calculated and marked as A_(i) .

Step 2.2.1.3.3.3: if A_(i) =|J_(i)|, S_(i)* is assigned to S_(i)* , and C_(i)* is assigned to C_(i)*. Thereafter, step 2.2.1.3.3.1 is carried out; otherwise, step 2.2.1.3.4 is carried out.

Step 2.2.1.3.4: the selected coverage opportunity set S_(i)* and the coverage mode set C_(i)* are outputted. The maximum value of the end time of the coverage opportunity in S_(i)* is used as the completion time F_(i).

The first phase of the two-phase heuristic algorithm in step 2.2.1.3 utilizes the rule based on dynamic greed. The coverage mode covering most valid grids is selected at each time. Once the coverage opportunity selects one coverage mode, the coverage opportunity may exit the subsequent selection activities. Once the grid is covered in the previous selection, the grid exits the subsequent selection activities. The method is simple, fast and effective and very suitable for rapidly obtaining the desirable feasible coverage solution under the current environment. At the second phase, on the basis of the solution obtained at the first phase, after the selected coverage opportunity with maximum end time is deleted, the complete coverage solution is re-planned to search a better solution space of the solution.

Step 2.2.2: the coverage opportunity sets {S_(i)*|i=1, 2, . . . , |N|} selected by all subregions are combined to obtain an overall selected coverage opportunity set S*. The coverage solutions {C_(i)*|i=1, 2, . . . , |N|} of all subregions are combined to obtain an overall coverage solution C* of the to-be-observed rectangular region R. The maximum value of the completion time {F_(i)|i=1, 2, . . . , |N|} corresponding to the coverage solutions of all subregions is used as the overall completion time F of the to-be-observed rectangular region R.

The solution space of the allocation solution is searched below by utilizing the simulated annealing method to obtain a better allocation solution, so that the completion time corresponding to the formed coverage solution is minimized, which specifically includes step 2.3-step 2.9:

Step 2.3: maximum iteration times k and initial annealing temperature T₁ are set, and the iteration times of the current outer cycle is set as k=1.

Step 2.4: if k>k, step 2.9 is carried out; otherwise, step 2.5 is carried out.

Step 2.5: the insert operator and the exchange operator are used to construct a plurality of neighbor allocation solutions of the current allocation solution Y.

Step 2.5.1: the number of the to-be-constructed neighbor allocation solutions is set as m, and the number of the constructed neighbor allocation solutions is set as m, and m=0.

Step 2.5.2: the insert operator is used to construct the neighbor allocation solution.

Step 2.5.2.1: an idle coverage opportunity ŝ is found from the unselected coverage opportunity set S′(Y)={S/S*(Y)} under the current allocation solution Y. The end time β(ŝ) of the idle coverage opportunity ŝ is smaller than the overall completion time F(Y) but greater than the completion time F_(y) _(ŝ) (Y) that the idle coverage opportunity ŝ is allocated to the subregion y_(ŝ). If the idle coverage opportunity cannot be found, step 2.5.2.4 is carried out; otherwise, step 2.5.2.2 is carried out.

Step 2.5.2.2: m+1 is assigned to m, and the idle coverage opportunity ŝ in the current allocation solutions Y is randomly allocated to one subregion, and the end time β(ŝ) of the idle coverage opportunity ŝ is smaller than the completion time of the allocated subregion, thereby forming a neighbor allocation solution Y _(m).

Step 2.5.2.3: the idle coverage opportunity ŝ is removed from the coverage opportunity set S′(Y), and then step 2.5.2.1 is carried out.

Step 2.5.2.4: all neighbor allocation solutions obtained by utilizing the insert operator are outputted.

Step 2.5.3: the exchange operator is used to construct the neighbor allocation solution.

Step 2.5.3.1: if m≥m, step 2.5.3.3 is carried out; otherwise, step 2.5.3.2 is carried out.

Step 2.5.3.2: m+1 is assigned to m. Two subregions n_(i) and n_(j) are randomly selected from the subregion set N. One coverage opportunity s_(i) is randomly selected from the coverage opportunity set {s|y_(s)=n_(i)} allocated to the subregion n_(i). One coverage opportunity s_(j) is randomly selected from the coverage opportunity set {s|y_(s)=n_(j)} allocated to the subregion n_(j). s_(i) and s_(j) are exchanged in the current allocation solution Y, that is, n_(j) is assigned to y_(s) _(i) , and n_(i) is assigned to y_(s) _(j) , thereby forming a neighbor coverage solution Y _(m), and then step 2.5.3.1 is carried out.

Step 2.5.3.3: all neighbor allocation solutions obtained by using the exchange operator are outputted.

Step 2.6: the coverage opportunity set selected under each neighbor allocation solution, desirable coverage mode set and corresponding completion time are calculated separately, and the neighbor allocation solution with minimum completion time is selected as the optimal neighbor allocation solution Y*. The specific steps refer to step 2.2 and are not repeated here.

Step 2.7: The completion time of the optimal neighbor allocation solution Y* is compared with the completion time of the current allocation solution Y. If F(Y*)≤F(Y) the optimal neighbor allocation solution Y* is assigned to the current allocation solution Y, and then step 2.8 is directly carried out; otherwise, step 2.8 is carried out after the step 2.7.1-step 2.7.2.

Step 2.7.1: an inferior solution accept rate

$P_{k} = {\exp \left\{ \frac{{F(Y)} - {F\left( {\overset{\_}{Y}}^{*} \right)}}{T_{k}} \right\}}$

at k^(th) iteration is calculated.

Step 2.7.2: a random number r_(k) of 0-1 is generated at k^(th) times. If r_(k)≤p_(k), the optimal neighbor allocation solution Y* is assigned to the current allocation solution Y.

Step 2.8: k+1 is assigned to k, and the annealing temperature T_(k)=T_(k-1)×λ is updated. λ is a cooling coefficient and is a fixed value less than 1 and greater than 0, and then step 2.4 is carried out.

The annealing temperature represents randomness for searching the solution space. If the requirement for the final solution quality is high, λ may be set to be large, and the iteration times may be increased. If the rapid convergence is required, λ is set to be small.

Step 2.9: the current allocation solution Y as well as the overall selected coverage opportunity set S*(Y) under the allocation solution, better coverage modes C*(Y) of each coverage opportunity and corresponding completion time F(Y) are outputted. 

We claim:
 1. A partition satellite mission planning method for earliest completion of regional target coverage, wherein the partition satellite mission planning method comprises: step 1: partitioning a to-be-observed rectangular region R; step 2: allocating observation resources to different regions, and selecting a certain coverage opportunity and a coverage mode thereof to minimize the completion time of a formed coverage solution, comprising: step 2.1: randomly allocating the coverage opportunities to each subregion to form an initial current allocation solution Y={y_(s)|∀s∈S}, wherein the value of y_(s) represents the subregion where the coverage opportunity s is allocated, and S represents a coverage opportunity set; step 2.2: calculating the coverage opportunity set S*(Y) selected under the current allocation solution Y, a better coverage mode set C*(Y) and corresponding completion time F(Y); step 2.3: setting maximum iteration times k and initial annealing temperature T₁, and setting the iteration times of the current outer cycle as k=1; step 2.4: if k>k, carrying out step 2.9; otherwise, carrying out step 2.5; step 2.5: using an insert operator and an exchange operator to construct a plurality of neighbor allocation solutions of the current allocation solution Y; step 2.6: calculating the coverage opportunity set selected under each neighbor allocation solution, better coverage mode set and corresponding completion time separately, and selecting the neighbor allocation solution with minimum completion time as the optimal neighbor allocation solution Y*; step 2.7: comparing the completion time of the optimal neighbor allocation solution Y* with the completion time of the current allocation solution Y; if F(Y*)≤F(Y), assigning the optimal neighbor allocation solution Y* to the current allocation solution Y, and then carrying out step 2.8 directly; otherwise, carrying out step 2.8 after the step 2.7.1-step 2.7.2; step 2.7.1: calculating an inferior solution accept rate $P_{k} = {\exp \left\{ \frac{{F(Y)} - {F\left( {\overset{\_}{Y}}^{*} \right)}}{T_{k}} \right\}}$ at k^(th) iteration; step 2.7.2: generating a random number r_(k) of 0-1 at k^(th) times; if r_(k)≤p_(k), assigning the optimal neighbor allocation solution Y* to the current allocation solution Y; step 2.8: assigning k+1 to k, and updating the annealing temperature T_(k)=T_(k-1)×λ, wherein λ is a cooling coefficient and is a fixed value less than 1 and greater than 0; and then carrying out step 2.4; step 2.9: outputting the current allocation solution Y as well as the overall selected coverage opportunity set S*(Y) under the allocation solution, better coverage modes c*(Y) of each coverage opportunity and corresponding completion time F(Y).
 2. The partition satellite mission planning method according to claim 1, wherein the step 1 is carried out as follows: step 1.1: setting length L_(p) and width W_(p) of each subregion; step 1.2: calculating the length L_(R) and width W_(R) of the to-be-observed rectangular region R according to the coordinates of four vertexes of the to-be-observed rectangular region R; step 1.3: dividing L_(p) by L_(R) to obtain length section number D_(L) and remainder R_(L) of the subregion; if the remainder ${R_{L} \geq {\frac{D_{L}}{{2D_{L}} + 1}L_{p}}},$ assigning D_(L)+1 to D_(L); dividing W_(p) by W_(R) to obtain width section number D_(W) and remainder R_(W) of the subregion, and if the remainder ${R_{W} \geq {\frac{D_{W}}{{2D_{W}} + 1}W_{p}}},$ assigning D_(W)+1 to D_(W); step 1.4: dividing L_(R) by D_(L) to obtain length L_(f) of the final subregion, and dividings W_(R) by D_(W) to obtain width W_(f) of the final subregion; step 1.5: according to the length L_(f) and width W_(f) of the final subregion, equally partitioning the to-be-observed rectangular region R, thereby obtaining the coordinates of four vertexes of each subregion, comprising: step 1.5.1: if i=1, using the left upper corner coordinate LU(R) of the to-be-observed rectangular region R as the left upper corner coordinate LU(n_(i)) of the i^(th) sub regionn n_(i); step 1.5.2: adding an abscissa value of the left upper corner coordinate LU(n_(i)) of the i^(th) subregion n_(i) and L_(f) to obtain the right upper corner coordinate RU(n_(i)) thereof; adding an ordinate value of the left upper corner coordinate LU(n_(i)) and W_(f) to obtain the left lower corner coordinate LD(n_(i)) thereof; and adding the ordinate value of the right upper corner coordinate RU(n_(i)) and W_(f) to obtain the right lower coordinate RD(n_(i)); step 1.5.3: if the right lower corner coordinate RD(n_(i)) of the i^(th) subregion n_(i) is equal to the right lower corner coordinate RD(R) of the to-be-observed rectangular region R, carrying out step 1.5.6; otherwise, carrying out step 1.5.4; step 1.5.4: if the abscissa value of the right upper corner coordinate RU(n_(i)) of the i^(th) subregion n_(i) is equal to the abscissa value of the right upper corner coordinate RU(R) of the to-be-observed rectangular region R, using LD(n_((i-D) _(L) ₊₁₎) as the left upper corner coordinate LU(n_(i+1)) of the (i+1)^(th) subregion n_(i+1); otherwise, using the right upper corner coordinate RU(n_(i)) of the i^(th) subregion n_(i) as the left upper corner coordinate LU(n_(i+1)) of the (i+1)^(th) subregion n_(i+1); step 1.5.5: after i+1 is assigned to i, carrying out step 1.5.2 to construct the coordinates of four vertexes of a next subregion; step 1.5.6: ending the partitioning, and outputting the coordinates of four vertexes of each subregion; step 1.6: obtaining a set N of the subregions, and the total number of the subregions of |N|=D_(L)×D_(W).
 3. The partition satellite mission planning method according to claim 1, wherein the step 2.5 is carried out as follows: step 2.5.1: setting the number of the to-be-constructed neighbor allocation solutions as m, and setting the number of the constructed neighbor allocation solutions as m, and m=0; step 2.5.2: using the insert operator to construct the neighbor allocation solution, comprising step 2.5.2.1-step 2.5.2.4: step 2.5.2.1: finding an idle coverage opportunity ŝ from the unselected coverage opportunity set S′(Y)={S/S*(Y)} under the current allocation solution Y, wherein the end time β(ŝ) of the idle coverage opportunity ŝ is smaller than the overall completion time F(Y) but greater than the completion time F_(y) _(ŝ) (Y) that the idle coverage opportunity ŝ is allocated to the subregion y_(ŝ); if the idle coverage opportunity cannot be found, carrying out step 2.5.2.4; otherwise, carrying out step 2.5.2.2; step 2.5.2.2: assigning m+1 to m, and randomly allocating the idle coverage opportunity ŝ in the current allocation solutions Y to one subregion, wherein the end time β(ŝ) of the idle coverage opportunity ŝ is smaller than the completion time of the allocated subregion, thereby forming a neighbor allocation solution Y _(m); step 2.5.2.3: removing the idle coverage opportunity ŝ from the coverage opportunity set S′(Y), and then carrying out step 2.5.2.1; step 2.5.2.4: outputting all neighbor allocation solutions obtained by utilizing the insert operator; step 2.5.3: using the exchange operator to construct the neighbor allocation solution, comprising step 2.5.3.1-step 2.5.3.3: step 2.5.3.1: if m≥m, carrying out step 2.5.3.3; otherwise, carrying out step 2.5.3.2; step 2.5.3.2: assigning m+1 to m; randomly selecting two subregions n_(i) and n_(j) from the subregion set N; randomly selecting one coverage opportunity s_(i) from the coverage opportunity set {s|y_(s) n_(i)} allocated to the subregion n_(i); randomly selecting one coverage opportunity s_(j) from the coverage opportunity set {s|y_(s)=n_(j)} allocated to the subregion n_(j); exchanging s_(i) and s_(j) in the current allocation solution Y, that is, assigning n_(j) to y_(s) _(i) , and assigning n_(i) to y_(s) _(j) , thereby forming a neighbor coverage solution Y _(m), and then carrying out step 2.5.3.1; step 2.5.3.3: outputting all neighbor allocation solutions obtained by using the exchange operator.
 4. The partition satellite mission planning method according to claim 1, wherein the step 2.2 and the step 2.6 of calculating the selected coverage opportunity set under the given allocation solution, a better coverage mode set selected for each coverage opportunity and corresponding completion time are carried out as follows: step a: using different subregions and the set of coverage opportunities allocated to different subregions as inputs to construct grids, and using a two-phase heuristic algorithm to calculate the coverage opportunity set selected for local part of different subregions as well as the better coverage solution and corresponding completion time; step a.1: dividing the i^(th) subregion n_(i) into a plurality of consistent square grids to obtain a grid set J_(i); step a.2: according to the divided grid set J_(i), generating a longest basic coverage mode for the coverage opportunity set S_(i) allocated to the i^(th) subregion n_(i), and forming a main set C_(i) of the optional coverage modes, wherein the longest basic coverage mode refers to that the length is equal to the maximum length r_(s), and the following two types of grids exist: {circle around (1)} one type of grids is completely covered by the coverage mode, that is, the four vertexes are within the area covered by the coverage mode; {circle around (2)} one grid has a vertex located on a left border of the coverage mode and is named left grid, and the other grid has a vertex located on an upper border of the coverage mode and is named upper grid, wherein the left grid and the upper grid can be overlapped; step a.3: by utilizing the grid set J_(i) separated from i^(th) subregion n_(i) as well as the coverage opportunity set S_(i)={s|y_(s)=n_(i), s∈S} allocated to i^(th) subregion n_(i) and optional coverage mode set C_(i) as inputs, determining the coverage opportunity subset S_(i)* selected by the i^(th) subregion n_(i) by the two-phase heuristic algorithm, and selecting the coverage mode of each coverage opportunity, so that the better coverage mode set C_(i)* is obtained, and the completion time F_(i) thereof can be obtained at the same time; step b: combining the coverage opportunity sets {S_(i)*|i=1, 2, . . . , |N|} selected by all subregions to obtain an overall selected coverage opportunity set S*; combining the coverage solutions {C_(i)*i=1, 2, . . . , |N|} of all subregions to obtain an overall coverage solution C* of the to-be-observed rectangular region R; and using the maximum value of the completion time {F_(i)|i=1, 2, . . . , |N|} corresponding to the coverage solutions of all subregions as the overall completion time F of the to-be-observed rectangular region R.
 5. The partition satellite mission planning method according to claim 4, wherein the step a.2 of generating a longest basic coverage mode for the coverage opportunity set according to the grid set and forming a main set of the optional coverage modes is carried out as follows: step a.2.1: if S_(i)=Ø, carrying out step a.2.4; otherwise, carrying out step a.2.2; step a.2.2: selecting one of coverage opportunities s∈S_(i) and according to the type of the coverage opportunity s, using methods in step a.2.2.1 to step a.2.2.10 to construct a longest basic coverage mode set C_(s) based on the grid set J_(i) for the coverage opportunity s; step a.2.2.1: judging the type of the coverage opportunity s according to the sub-satellite point track o_(s) of the coverage opportunity s, and the types including a right lower corner oblique type and a left lower corner oblique type; step a.2.2.2: making the set as C_(s)=Ø; step a.2.2.3: selecting subsets J_(s) from the grid set J_(i), and using each grid in the subsets J_(s) as the left grid of the longest basic coverage mode of the coverage opportunity s; step a.2.2.3.1: setting J_(s)=Ø; step a.2.2.3.2: traversing the set J_(i), and calculating the camera roll angle ${\sigma_{s}(u)} = {{\arctan \frac{d_{s}^{u}}{h_{s}}} - \frac{w_{s}}{2}}$ of the coverage opportunity s when each grid u is used as the left grid of the coverage mode, wherein if the coverage opportunity is the left lower corner oblique type, d_(s) ^(u) is a distance between the left lower corner vertex LD(u) of the grid u and the sub-satellite point track o_(s) of the coverage opportunity s; if the coverage opportunity is the left lower corner oblique type, d_(s) ^(u) is a distance between the left upper corner vertex LU(u) of the grid u and the sub-satellite point track o_(s) of the coverage opportunity s; if |σ_(s) (u)≤v_(s), adding the grid u into the subset J_(s) of the left grid; step a.2.2.4: if J_(s)=Ø, carrying out step a.2.2.10; otherwise, after one grid p in J_(s) is selected, carrying out step 2.2.1.2.2.5; step a.2.2.5: according to the grid p, selecting the grid subset J_(s) ^(p) from the grid set J_(i), and when p is used as the left grid, using each grid in J_(s) ^(p) as the upper grid to form a longest basic coverage mode of the coverage opportunity s together with P; step a.2.2.6: if J_(s) ^(p)=Ø, carrying out step a.2.2.9; otherwise, after one grid q in J_(s) ^(p) is selected, carrying out step a.2.2.7; step a.2.2.7: based on the grids p and q, constructing the longest basic coverage c_(s) (p,q) belonging to the coverage opportunity s; step a.2.2.8: adding c_(s)(p,q) into C_(s), removing the grid q from J_(s) ^(p), and then carrying out step a.2.2.6; step a.2.2.9: removing the grid P from J_(s), and then carrying out step a.2.2.4; step a.2.2.10: outputting the coverage mode set C_(s); step a.2.3: removing the coverage mode opportunity s from S_(i), and then carrying out step a.2.1; step a.2.4: if C_(i)={c_(s)|✓c_(s)∈C_(s), ∀s∈S_(i)}, outputting the coverage mode set C_(i).
 6. The partition satellite mission planning method according to claim 5, wherein the step a.2.2.5 is carried out as follows: step a.2.2.5.1: setting J_(s) ^(p)=Ø; step a.2.2.5.2: calculating a strip width η_(s)(p) of the coverage mode when the grid p is used as the left grid by using the formula (1): $\begin{matrix} {{\eta_{s}(p)} = \left\{ \begin{matrix} {{d_{s}^{p} - {h_{s} \cdot {\tan \left( {{\arctan \frac{d_{s}^{p}}{h_{s}}} - w_{s}} \right)}}},} & {{\arctan \frac{d_{s}^{p}}{h_{s}}} > w_{s}} \\ {{d_{s}^{p} + {h_{s} \cdot {\tan \left( {w_{s} - {\arctan \frac{d_{s}^{p}}{h_{s}}} - w_{s}} \right)}}},} & {{\arctan \frac{d_{s}^{p}}{h_{s}}} \leq w_{s}} \end{matrix} \right.} & (1) \end{matrix}$ step a.2.2.5.3: traversing the set J_(i); checking whether one grid t satisfies the following four conditions; if all conditions are satisfied, adding the grid t into J_(s) ^(p); otherwise, checking the next grid; if the coverage opportunity is the right lower corner oblique type, the four conditions are: {circle around (1)} making the straight line L_(p)∥o_(s) through the left lower corner vertex of the grid p, making the straight line L_(t)∥o_(s) through the right upper corner vertex of t, and ensuring that the distance between the straight lines L_(p) and L_(t)′ do not exceed the width η_(s)(p) of the strip; {circle around (2)} making the straight lines L_(p)′⊥o_(s) and L_(t) ⊥o_(s) respectively through the left upper corner vertexes of the grids p and t, and enabling L_(t) to be above L_(p)′; {circle around (3)} making the straight line L _(t) to be through the left lower corner vertex of the grid t, and enabling the straight line L _(t) to be at the right side of L_(p); {circle around (4)} making the straight line L _(p)⊥o_(s) through the right lower corner vertex of the grid p, and ensuring that the distance between the straight lines L _(p) and L_(t) do not exceed the length l_(s) of the strip, wherein “∥” indicates a parallel relation, and “⊥” indicates a vertical relation; if the coverage opportunity is the left lower corner oblique type, the four conditions are: {circle around (1)} making the straight line L_(p)∥o_(s) through the left upper corner vertex of the grid p, making the straight line L_(t)′∥o_(s) through the right upper corner vertex of t, and ensuring that the distance between the straight lines L_(p) and L_(t)′ do not exceed the width η_(s)(p) of the strip; {circle around (2)} making the straight lines L_(p)′⊥o_(s) and L_(t)⊥o_(s) respectively through the right upper corner vertexes of the grids p and t, and enabling L_(t) to be above L_(p)′; {circle around (3)} making the straight line L _(t)∥o_(s) through the left upper corner vertex of the grid t, and enabling the straight line L _(t) to be at the right side of L_(p); 4 making the straight line L _(p)⊥o_(s) through the left lower corner vertex of the grid p, and ensuring that the distance between the straight lines L _(p) and L_(t) do not exceed the length l_(s) of the strip, wherein “∥” indicates a parallel relation, and “⊥” indicates a vertical relation.
 7. The partition satellite mission planning method according to claim 5, wherein the step a.2.2.7 is carried out as follows: step a.2.2.7.1: if the coverage opportunity is the right lower corner oblique type, making the straight line L_(p)∥o_(s) through the left lower corner vertex LD(p) of the grid p, and making the straight line L_(q)⊥o_(s) through the left upper corner vertex LU(q) of the grid q; if the coverage opportunity is the left lower corner oblique type, making the straight line L_(p)∥o_(s) is made through the left upper corner vertex LU(p) of the grid p, and making the straight line L_(q)⊥o_(s) through the right upper corner vertex RU(q) of the grid q; step a.2.2.7.2: calculating the strip width η_(s)(p) of the coverage mode when the grid p is used as the left grid; step a.2.2.7.3: making the straight line L _(p) so as to satisfy: {circle around (1)} L _(p)∥L_(p); {circle around (2)} the distance between L _(p) and L_(p) is equal to η_(s)(p); {circle around (3)} L _(p) is at the right side of L_(p); step a.2.2.7.4: making the straight line L _(q) so as to satisfy: {circle around (1)} L _(q)∥L_(q); {circle around (2)} the distance between L _(q) and L_(q) is equal to l_(s); L _(q) is at the lower side of L_(q); step a.2.2.7.5: using four straight lines of L_(p), L _(p), L_(q), L _(q) as four sides of the rectangle, and using intersections of the four straight lines as four vertexes of the rectangle to construct the coverage mode c_(s)(p,q).
 8. The partition satellite mission planning method according to claim 4, wherein the step a.3 is carried out as follows: step a.3.1: by utilizing the grid set J_(i), the coverage opportunity set S_(i) and the corresponding coverage mode set C_(i) as inputs, calculating the selected coverage opportunity subsets S_(i)*, the selected coverage mode sets C_(i)* and the number A_(i) of corresponding coverable grids with the heuristic algorithm; step a.3.2: if A_(i)<|J_(i)|, with |J_(i)| indicating the number of elements in the set J_(i), carrying out step a.3.4; otherwise, carrying out step a.3.3; step a.3.3: optimizing S_(i)*, and C_(i)*, by using a promotion algorithm; step a.3.3.1: removing one coverage opportunity $s_{m} = {\arg {\max\limits_{\forall{s \in S_{i}^{*}}}{\beta (s)}}}$ with maximum end time in S_(i)* from S_(i)* to obtain the coverage opportunity set S_(i) ; step a.3.3.2: by utilizing the grid set J_(i), the coverage opportunity set S_(i) and the optional coverage mode set C_(i) ={c|∀c∈C_(s), ∀s∈S_(i) } corresponding to S_(i) as inputs, calculating the selected coverage opportunity subsets S_(i)*, the selected coverage mode set C_(i)*, and the number A_(i) of the corresponding coverable grids with the heuristic algorithm; step a.3.3.3: if A_(i) =|J_(i)|, assigning S_(i)* to S_(i)*, and assigning C_(i)* to C_(i)*; carrying out step a.3.3.1; otherwise, carrying out step a.3.4; step a.3.4: outputting the selected coverage opportunity set S_(i)* and the coverage mode set C_(i)*; and using the maximum value of the end time of the coverage opportunity in S_(i)* as the completion time F_(i).
 9. The partition satellite mission planning method according to claim 8, wherein the step a.3.1 is carried out as follows: step a.3.1.1: sorting all coverage opportunities in S_(i) in a non-descending order according to the end time to obtain; ${\overset{\rightharpoonup}{S_{i}} = \left\{ {s_{1},s_{2},\ldots,s_{|\overset{\rightharpoonup}{s_{i}|}}}\  \right\}};$ step a.3.1.2: setting the state of all grids in J_(i) as “uncovered”; g=0, wherein g indicates the subscript of the g^(th) coverage mode in the coverage mode set C_(i)* ; C_(i)*=Ø; step a.3.1.3: assigning g+1 to g; extracting the coverage opportunity s_(g) from {right arrow over (S_(i))}; traversing the optional coverage mode set C_(s) _(g) ; for each coverage mode c_(s) _(g) , selecting the grid capable of being completely covered by c_(s) _(g) and with the state of “uncovered” from the grid set J_(i) and marking as valid grid sets V(c_(s) _(g) ), and marking the number as |V(c_(s) _(g) )|; step a.3.1.4: selecting the coverage mode $c_{s_{g}}^{*} = \left. {\arg \max\limits_{\forall{c_{s_{g}} \in C_{s_{g}}}}} \middle| {V\left( c_{s_{g}} \right)} \right|$ covering the most valid grids from C_(s) _(g) to serve as the coverage mode selected by the coverage opportunity s_(g); after c_(s) _(g) * is added into C_(i)*, changing the state of all grids in the valid grid sets V(c_(s) _(g) ) of the coverage mode c_(s) _(g) * as “covered”; step a.3.1.5: calculating the number of all grids with the state of “covered” in J_(i) and marking as A_(i); step a.3.1.6: if A_(i)=|J_(i)| or g=|{right arrow over (S_(i))}|, carrying out step a.3.1.7; otherwise, carrying out step a.3.1.3; step a.3.1.7: outputting the coverage opportunity subset S_(i)*={s₁, s₂, . . . , s_(g)} as well as the coverage mode set C_(i)* thereof and the number A_(i) of covered grids.
 10. The partition satellite mission planning method according to claim 8, wherein the step a.3.3.2 is carried out as follows: step a.3.3.2.1: setting the state of all grids in J_(i) as “uncovered”, and C_(i)*=Ø; step a.3.3.2.2: g=1, wherein g indicates the subscript of the g ^(th) coverage mode in the coverage mode set C_(i)*; step a.3.3.2.3: if g>|C_(i) |, carrying out step a.3.3.2.5; otherwise, carrying out step a.3.3.2.4; step a.3.3.2.4: selecting the g ^(th) coverage mode in C_(i) , and selecting the grids capable of being completely covered by the g ^(th) coverage mode and with the state of “uncovered” from the grid set J_(i) and marking as valid grid set V(g), and marking the number as |V(g)|; assigning g+1 to g; and carrying out step a.3.3.2.3; step a.3.3.2.5: if C_(i) =Ø, carrying out step a.3.3.2.8; otherwise, carrying out step a.3.3.2.6; step a.3.3.2.6: selecting the coverage mode ${\overset{\_}{g}}^{*} = \left. {\arg \max\limits_{\forall{\overset{\_}{g} \in \overset{\_}{C_{i}}}}} \middle| {V\left( \overset{\_}{g} \right)} \right|$ with most valid grids in C_(i) to serve as the coverage mode selected by the corresponding coverage opportunity s*; adding g* into C_(i)*, and changing the state of all grids in the valid grid set V(g*) of the coverage mode g* to be “covered”; step a.3.3.2.7: removing the coverage mode set C _(s*) corresponding to the coverage opportunity s* from the optional coverage mode set C_(i) , and then carrying out step a.3.3.2.2; step a.3.3.2.8: outputting the coverage mode set C_(i)* selected by each coverage opportunity, i.e. the final coverage solution; and calculating the number of all “covered” grids in J_(i) and marking as A_(i) . 